Physics of the liquid-hquid critical point 
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Within the inherent structure {IS) thermodynamic formalism introduced by Stillinger and Weber 
[F. H. Stillinger and T. A. Weber, Phys. Rev. A 25, 978 (1982)] we address the basic question of the 
physics of the liquid-liquid transition and of density maxima observed in some complex liquids such 
as water by identifying, for the first time, the statistical properties of the potential energy landscape 
(PEL) responsible for these anomalies. We also provide evidence of the connection between density 
anomalies and the liquid-liquid critical point. Within the simple (and physically transparent) model 
discussed, density anomalies do imply the existence of a liquid-liquid transition. 
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Water, the most important liquid for life, belongs to 
a class of liquids for which the isobaric temperature de- 
pendence of the density has a maximum. In contrast to 
ordinary liquids, water expands on cooling below AC at 
atmospheric pressure. This density anomaly is associated 
with other thermodynamic anomalies, such as minima in 
the compressibility along isobars and an increase of the 
specific heat on cooling 1]. Recent fascinating studies 
have attempted to connect these anomalies to the pres- 
ence of two different liquid structures, separated at low 
temperature by a line of first order transitions, ending 
in a second order critical point H, |^ Q . In the case of 
water, this novel critical point would be located in an 
experimentally inaccessible region [^. Despite this un- 
favorable location, the postulated presence of this crit- 
ical point provides a framework for interpreting |Q| not 
only features of the liquid state but also the phenomenon 
of polyamorphism and the first-order like transition be- 
tween polyamorphs 0- In this Letter we aim at iden- 
tifying the statistical properties of the potential energy 
landscape (PEL) responsible for the density maxima and 
the connection to the physics of the liquid-liquid transi- 
tion. 

The study of the statistical properties of the PEL — 
i.e. of the number, shape and depth of the basins compos- 
ing the potential energy hypersurface — is under tremen- 
dous development. Theoretical approaches, based on the 
seminal work of Stillinger and Weber combined ei- 
ther with calorimctric data ^ .IO lJ I J or with analysis of 
extensive numerical simulations |l2l lisl are starting 
to provide precise estimates of the statistical properties 
of the PEL in several materials and models for liquids. 
A simple model for the statistical properties of the land- 
scape, supported by recent numerical studies p^Hslll^ . 
can be built on the basis of the two following hypotheses: 

1. a g aussian distribution for Q.{eis)deis jl^L \v\ \\A 
lig . |2flj, the number of distinct basins of energy 
depth e/s between ejs and ejs + deis, i.e. 



VL{eis)deis = 



V2 



-de IS 



(1) 



Here e"^ counts the total number of PEL basins, 
{N being the number of molecules) Eg is the char- 
acteristic energy scale and cr^ is the variance of the 
distribution. A gaussian distribution is suggested 
by the central limit theorem, since in the absence of 
diverging correlation, ejs can be written as the sum 
of the IS energy of several independent subsystems. 

A (multi dimensional) parabolic shape of the PEL 
basins. This hypothesis is equivalent to assuming 
that the T-dependence of basin exploration is gov- 
erned by an harmonic vibrational free energy fyu,. 
For a system of rigid water molecules 



6N~3 



fvib{eis,T) ^ ksT ^ ln{phuji{eis))] 



(2) 



where ?i is Planck's constant, /3 = l/{kBT) and 
ks the Boltzmann constant and uJi is the fre- 
quency of the z-th normal mode. As suggested by 
numerical evidence 0, 0, 0|, the relation be- 
tween the parabolic basin volume and the basin 
depth is assumed to be linear and the linearity 
is quantified by two parameters a and b by writ- 
ing Yl^i=i^ '^^{'^ii^is)) = a, + beis- More accurate 
descriptions for fvibi^is^T), incorporating anhar- 
monic corrections, can be employed '2 1] if accuracy 
in the model is required. In this work, which fo- 
cuses on the connection between landscape features 
and thermodynamic scenarios, the simple harmonic 
expression is considered. Accounting for anhar- 
monic corrections does not affect the main conclu- 
sion of this work. 

For this harmonic gaussian landscape model, the three 
parameters a, Eq, (modeling the statistical properties 



2 



of the PEL) and the two parameters a and h (modehng 
the relation between basin shape and depth) fully deter- 
mine the free energy of the liquid. The corresponding 
equation of state (EOS), expressed in term of the V- 
derivative of the landscape properties, is 

P(r, V) = V const + TV^ + T-'v^^^ , (3) 

where VconstiV) - -^[^o - ba% V^{V) = R^[a ~ 
a-bEo + and V^^^ {V) = ^[a^ /2R]. 

Along an isochore, the high T behavior is fixed by the 
linear term TV , while the low T behavior is controlled 
by the T^^ term T^^V^ ^ . From Eq. |21 one can also 
conclude that, in the harmonic gaussian landscape, P 
along an isochore either is monotonically increasing with 
T (if V^^^ < 0, simple liquid cases), or it has a minimum 

at T = ^'P^^^{y)/'P^{y) (if > 0, liquids with 

density anomalies) 23] . 

Eq. O offers the possibility of understanding the land- 
scape parameters responsible for density anomalies. In- 
deed, a Maxwell relation states that a density maximum 
state point (i.e. a point where dV/dT\p = 0) is simulta- 
neously a point at which the T-dependence of P along an 
isochore has a minimum {dP/dT\v = 0). Hence, the con- 
dition for the existence of density maxima, from Eq. |21 is 
that P^^^ > which, in the PEL formalism, corresponds 

to da^/dV > 0. Thus, in normal liquids, in the y-range 
where liquid states exist, da^/dV must always be neg- 
ative. On the contrary, liquids with density anomalies 
must be characterized by a V^-range where increases 
with V HI. 

The harmonic gaussian landscape also elucidates the 
relation between density anomalies and the existence of 
a liquid-liquid critical point. Fig. ^ schematically com- 
pares P{T, V) along three different isochores: at the vol- 
ume at which is a minimum, V — V„2 , and at two 
other ^-values, respectively below and above V^2 . In 
the harmonic gaussian landscape P along the V„2 iso- 
chore increases linearly with T, since, by construction, 
V^^^ {V) = and > 0. For V < V„2 , P decreases 
on cooling, since P^^_^ {V) < 0. The opposite trend is 
observed for V > V^2 . An isothermal cut of the three 
isochores shows the corresponding P{V) curves. At high 
T, where the linear T term is dominant, P{V) is mono- 
tonically decreasing with V, while at low T the P{V) 
isotherm is not monotonic, indicating the presence of a 
region of instability (negative compressibility). As a re- 
sult, one must conclude that a critical point exists at an 
intermediate T, which can be accessed if no other mech- 
anism (as for example crystallization or glass transition) 
preempts its observation. 

We now compare the harmonic gaussian landscape pre- 
dictions with an extensive study of the PEL of one of 




FIG. 1: Schematic explanation of the connection between a 
minimum in cr^(y) and a hquid-liquid critical point. Square, 
circle and diamond indicate respectively a V smaller than, 
equal to and larger than V„2^^^ . The central panel shows 
the corresponding isochores. Note that, as shown in Eq|3 
P increases (decreases) on cooling when da^/dV > (< 0). 
The bottom panels report the corresponding high and low T 
isotherms. At low T an unstable Van der Waals loop appears. 



the most widely studied models for water, the extended 
simple-point charge (SPC/E) H^. SPC/E is able to 
reproduce water density anomalies. The PEL for this 
model has been studied in some detail previously [lRl2^ . 
We simulated a system of 216 water molecules with classi- 
cal molecular dynamics (MD) in the NVE ensemble. We 
studied 45 different state points in supercooled states, for 
times longer than several a-relaxation times. Long range 
forces have been modeled via the reaction field technique. 
The integration time step was 1.0 fs. Results are aver- 
aged over 100 different independent trajectories for each 
state point. From each of the 100 trajectories we ex- 
tracted about 30 configurations and calculated, via con- 
jugate gradient techniques with lO"""^^ tolerance, the lo- 
cal minima (the so-called inherent structures) and their 
energies e/5. The vibrational density of states in the 
local minima has been calculated by diagonalizing the 
Hessian (the matrix of second derivatives of the poten- 
tial). This procedure produces 3000 (30 x 100) distinct 
minima for each of the 45 studied state points, allow- 
ing the statistical error to be smaller than the symbol 
size. The resulting extremely accurate determination of 
Y^^=i ^ lii(wi(e/5)), of the T-dependence of (e/5) and of 
the configurational entropy Sconf |l5| — a measure of 
the logarithm of the number of explored PEL basins — 
is crucial when accurate ^-derivatives of the PEL param- 
eters are required, as in the present case. Procedures for 
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the evaluation of the PEL parameters for the harmonic 
gaussian model have been worked out 0, 0, 0, l2^ . 
From numerical evaluation of the T-dependence of (e/g) 
and Sconf it is possible (i) to provide evidence that the 
gaussian landscape correctly describes the SPC/E simu- 
lation data and (ii) calculate, with a fitting procedure, 
and Eq and a. The basin shape parameters a and b 
are calculated from the eis dependence of the vibrational 
density of states, evaluated at the IS, i.e. by a linear fit 
of ^ lii(wi(e7s)) vs e/g. The parameters Eq and 

are evaluated via a linear fit of (e/s(T)) vs T~^. The 
parameter a is calculated by fitting the T-dependence of 

SconfMM- 

Fig. 121 shows the T^-dependence of the landscape pa- 
rameters (a, Eq, (T^) for nine different volumes and 
contrasts it with the behavior characteristic of simple 
liquids |27|. The total number of states a decreases on 
compressing the system, a feature common to all sim- 
ple liquid models studied so far 0, . The V depen- 
dence of the energy scale Eq is also analogous to the one 
found in simple liquid models. Indeed Eq first decreases 
on compression, corresponding to the progressive sam- 
pling of the most attractive part of the potential, then 
it starts to increase due to the sampling of the repul- 
sive part of the potential at short intermolecular nearest 
neighbor distances. As expected on the basis of the pre- 
vious discussion, the ^-dependence of shows instead 
elements which are not observed in simple liquid models, 
where decreases monotonically on increasing V. In the 
SPC/E case, the variance shows a clear minimum around 
V^2 = 15 cm'^/mol (i.e. at density p = 1.2 g/cm'^) and 
hints of a maximum at F = 20 cm'^/mol (p = 0.9 g/cm^). 
Between 15 and 20 cm'^/mol, SPC/E water exhibits den- 
sity anomalies. The increase of cr^ for V > V^2 can be 
attributed to the fact that the system can explore new 
configurations, characterized by the presence of hydro- 
gen bonds. The formation of such bonds requires a large 
local volume and it is severely hampered at high density. 
These new accessible states widen the variance of the 
fl{eis) distribution, producing a range of 1^-values where 
da'^/dV is positive. When V has increased to about 20 
cm'^/mol the system has become mostly composed of a 
network of open linear hydrogen bonds and there are no 
more mechanisms available to increase the variance. 

The actual location of the critical point, which will 
depend on all landscape properties, can be calculated 
according to Eq. O For the SPC/E harmonic gaussian 
landscape, the resulting phase diagram is shown in Fig.O 
The density maxima locus retracing at low densities is 
consistent with the existence of a volume beyond which 
the liquid returns to normal (i.e. with da^/dV < 0). The 
T of maximum density (TMD) locus at large pressures 
terminates along the liquid-liquid spinodal, as predicted 
by thermodynamic consistency 01- Fig. 13 also show the 
TMD hne evaluated directly from the MD P{V, T) data. 
This curve, which properly includes anharmonic contri- 
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FIG. 2: Volume dependence of the PEL parameters for the 
SPC/E potential. Dashed lines indicates the trends expected 
for simple liquids 12^. 



butions, show that anharmonicities do not change the 
topology but merely shifts the TMD curve up about 25 
K in r and down 160 MPa in P. The location of the 
liquid- liquid critical point for the SPC/E potential is es- 
timated to be between w 140 and 175 K and between 
w 185 and 340 MPa, in agreement with previous predic- 
tions jH. 

One of the powerful features of the PEL formalism is 
the possibility of simultaneously evaluating the EOS and 
Sconf- This offers the unique possibility of evaluating, 
within the chosen landscape model, the theoretical limit 
of validity of the low T extrapolations. In the present 
case, the limit set by the condition Sconf = defines a 
Kauzmann curve Tk{P, V) below which the above EOS 
is no longer valid. The Tk{P,V) locus is also indicated 
in Fig. O The TMD locus crosses the Kauzmann lo- 
cus very close to its reentrance, as recently predicted by 
Speedy [sOl- It is interesting to observe that the loca- 
tion of the liquid-liquid critical point is within the glass 
state and hence technically does not exist in the har- 
monic gaussian approximation. When anharmonic con- 
tributions are taken into account, the critical point moves 
into the Sconf > region, suggesting that the liquid- 
liquid critical point could be a real feature of the SPC/E 
EOS. If this were the case, the possibility of observing 
a liquid-liquid critical point in this model would be only 
hampered by the extreme slowing down of the dynamics 
at low T and, in principle, by crystallization^lj. 

In this Letter we have shown that, in the harmonic 
gaussian landscape model, density anomalies are a suf- 
ficient, but not a necessary, condition for the existence 
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FIG. 3: Harmonic gaussian landscape phase diagram for the 
SPC/E potential, including the TMD line, the liquid-liquid 
critical point with the associated spinodal lines and the Kauz- 
mann locus. Lines below the Kauzmann locus have no phys- 
ical meaning, since they correspond to states with negative 
Sconf- They are drown here to visualize the termination of 
the TMD line and the relative position of the critical point 
respect to the Kauzmann locus. Squares indicate the TMD 
as evaluated directly from the MD calculations. As discussed 
in the text, adding the anharmonic contributions in the theo- 
retical calculation would move the liquid-liquid critical point 
into the physically meaningful Sconj > region. 

of a liquid- liquid critical point. Anharmonicities do not 
alter the picture and can be taken properly into account 
if necessary. Instead, the gaussian hypothesis is crucial, 
and models with a different distribution n{eis) may pro- 
duce density anomalies in systems without any liquid- 
liquid critical point [131 ■ Also, high-temperature fluid- 
fluid transitions|3,lj can not be described by this for- 
malism, which is applicable only to supercooled states. 
Still, the generality of the gaussian distribution, rooted 
in the central limit theorem strongly support the possi- 
bility that what we learn studying the gaussian landscape 
applies to water. 
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